clc;clear;close all;
m=0:.01:1;
[R,T]=ndgrid(m,-.6:.01:20);
x=1-(5/4*(1-mod(3.6*T, 2)).^2-1/4).^2/2;
P=pi*exp(-T/8)/2;
s=sin(P);
c=cos(P);
y=2*R.^2.*(1.3*R-1).^2.*s;
S=x.*(R.*s+y.*c);
X=S.*sin(T*pi);
Y=S.*cos(T*pi);
Z=R.*c-y.*s;

figure('color',[0,0,0])
subplot(2,2,1)
surf(X,Y,x.*Z,hypot(hypot(X,Y),Z*.9));
shading interp
axis equal off
m(2:3,:)=0;
colormap(m')

subplot(2,2,2)
surf(X,Y,x.*Z,hypot(hypot(X,Y),Z*.9));
shading interp
axis equal off
m(2:3,:)=0;
colormap(m')
view([90,0]);


subplot(2,2,3)
surf(X,Y,x.*Z,hypot(hypot(X,Y),Z*.9));
shading interp
axis equal off
m(2:3,:)=0;
colormap(m')
view([0,0]);


subplot(2,2,4)
surf(X,Y,x.*Z,hypot(hypot(X,Y),Z*.9));
shading interp
axis equal off
m(2:3,:)=0;
colormap(m')
view([0,90]);

